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VARIANCE REDUCTION BY IMPORTANCE SAMPLING AND THE METHOD OF 


SPLITTING IN MONTE CARLO CALCULATIONS 
by Burt M. Rosenbaum 
Lewis Research Center 

SUMMARY 

The two techniques of variance reduction that are considered are (1) importance 
sampling and (2) splitting and Russian roulette. Based on the value of the variance, op- 
timum biasing sampling procedures are investigated and it is determined when adjoint 
biasing yields the minimum variance. It is shown that the method of Russian roulette 
may lead to an increase, rather than a decrease, in variance. A short example illus- 
trates the methods used. 


INTRODUCTION 

Before the advent of the computer, when a problem involved a large number of 
members or participants (e. g. , when dealing with, say, a collection of molecules), 
analyses usually could not be carried out for the general case. Because of the complex- 
ities, analytical solutions could only be found in limiting situations or where simplifying 
assumptions could be made. In the regions where such assumptions could not be made, 
realistic theoretical analyses could not be accomplished and so-called educated guesses 
were resorted to. 

The computer enables an investigator to compare the theoretical behavior of his 
conceptual model with experimental data in the complicated intermediate region where 
standard analyses break down. When the Monte Carlo method is used, the possibilities 
of occurrence are encoded into the computer program and the behavior of a large number 
of sample particles is simulated by computer decisions as the sample particles are fol- 
lowed through the system. On the basis of the data thus generated average behavior 
patterns may be calculated. According to the information needed, the computer may be 
instructed to spew out local densities, total kinetic energy densities, heat and mass 
transfer rates, pressures, probabilities of penetration through a barrier, fission rates, 
chemical reaction rates, and so forth. 



In an attempt to decrease the computation time necessary for answers as well as to 
enable one to handle problems which originally would overload the computer capacity, 
techniques were evolved which established more efficient simulation processes than di- 
rect simulation; that is, better accuracy could be obtained for a given sample size. The 
systemization of such error- reducing procedures was due in large measure, to the work 
of H. Kahn (refs. 1 to 4). 

This report concerns itself with two of these techniques: (1) importance sampling 
and (2) splitting and Russian roulette. The theory has been extensively treated in the lit- 
erature and applications of the two techniques, used separately or in combination, 
abound in computer programs (refs. 5 to 13). In this report, both techniques are formu- 
lated in a unique manner and equations for the variance resulting from their use are de- 
rived. These equations are generalized to apply to any number of sets of random varia- 
bles and optimum procedures are developed. In addition, a brief example is posed and 
analyzed to illustrate the method. It is hoped that the formulation of the problem as pre- 
sented herein serves to clarify the basic concepts involved. 


IMPORTANCE SAMPLING 

The following problem is considered. Suppose we are given a function g dependent 

on three or more sets of random variables x = x., x 9 , . . . , x where x. represents all 
th l ci n i 

variables in the i set of random variables and we wish to determine the mean or ex- 
pectation value of g: 


E[g] = f ■ ■ ■ J g(x)f(x)dx 


( 1 ) 


In equation (1), the multivariate probability density function is denoted by f(5T) = f(x^, 

x 9 , . . . ,x ). When the relation 
l ’ n 

f(x) = f(x 1 ,x 2 )f(x 3 , . . ,,x n /x 1 ,x 2 ) (2) 


is used where the symbol f(x 3 , . . . . x /xj,x 2 ) is the conditional probability density 
function of the random variables x 3 , . . . ,x given that XpX 2 have taken on fixed val- 
ues, equation (1) can be written 


E[gl = If f(xp x 2 )dx 1 dx 2 I f g(x)f(x 3 , . . . ,x n /x 1 ,x 2 )dx 3 . . . dx n 

-If E[g/x 1 ,x 2 |f(x 1 ,x 2 )dx 1 dx 2 (3) 


2 


It can be seen, therefore, that Efg/x^x,], the expectation value of g given that 
and x 2 have been fixed, satisfies 

E[g/x 1 ,x 2 l = /•/ g(x)f(x 3 , . . . ,x n /x 1 ,x 2 )dx 3 - . ,dx n (4) 

The variance of g for fixed x^ and x 2 is given by 

4 /x V* 2 ‘f f (g(x) - Efg/x^x^j 2 f(x 3 ,. . . ,x n /x 1 ,x 2 )dx 3 . . .dx n (5) 

If we pick from a population distributed in accordance with the probability density 
function f(x^,x 2 ), then a weight function w(Xj,x 2 ) can be incorporated into the density 
function which may act to decrease the value of the variance (Xg/xpX 2 while keeping 
the value of E[g] invariant. The new probability density function f(Xj, x 2 )/w(x^, x 2 ) 
satisfies the relation 

f(x x ) 

1 dx 1 dx 9 = 1 (6) 

w(x r x 2 ) 

and equation (1) becomes 

f(x x . . . , x ) 

w(x r x 2 )g(x) dx 1 dx ? . . . dx (7) 

1 Z w( Xl ,x 2 ) 1 £ n 

The method just described whereby a weight function is employed is called "importance 
sampling" and, by equation (7), we see that the expectation value of g(x) when sampling 
from a population with probability density function f(x) is equal to the expectation value 
of w(xj,x 2 )g(x) when sampling from a population described by the density function 
f(x)/w(x r x 2 ). 

The variance of w(x^,x 2 )g(x) associated with the probability density function 
f(x )/w(Xp x 2 ) is given by 




3 



3 w(x 1 ,x 2 )g(x) 


f(x) 


w(x r x 2 ) 


J [ ( w<x l' x 2> g<x) ■ E fej) 2 *1 dx 2’ ■ -^n 

If w(Xj, x 2 )E^g 2 /x 1 , xJfCXj, x^dXj dx 2 - E 2 [g] 

E[w(x 1 ,x 2 )E[g 2 /x 1 ,x 2 


- E 2 [g] 


where 


g 2 (x)f(x 3 ,x 4 , . . . ,x n /x 1 ,x 2 )dx 3 - . ,dx n 

2 2 
= ff g /Xj,x 2 + E [g/x r x 2 ] 

The particular weight function w(Xj,x 2 ) that minimizes the variance of w(x^,x 2 )g(x) 
subject to the constraint given by equation (6) satisfies the equation 


( 9 ) 


E 


w(x x ,x 2 ) - __ b. 


ys^/xj, 


X' 


Y E [g 2 /x v x 2 \ 
and the minimum value of the variance is given by 


( 10 ) 


ft 

w(x 1 ,x 2 )g(x) 

f(x) 

II 

M 

tSD 

ft 

g 2 /x p x 2 


w(x 1 ,x 2 ) 





" E 2 [g] 


(ID 


This formulation, as mentioned before, is the method of importance sampling in the two 
variable sets x^ and x 2 . 

It may be stated here that if N measurements of w(xj,x 2 )g(x) were made, then the 
variance of the average is merely equal to the value given by equation (8) divided by N. 

The concept of stages is now introduced and it will be employed throughout this re- 
port. We assume that in the first stage we sample from the x^-distribution and in the 
second stage from the x 2 -distribution. We write 
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f(XpX 2 ) = f(x 1 )f(x 2 /x 1 ) 

(12) 

w(XpX 2 ) = uU^uCXpXg) 

(13) 

where u(x^) and u(Xp x 2 ) satisfy 


f f < x l) 

1 — Ldx = l 

J u(x i> 

(14a) 



f (x p/x i ) 
^-J-dX 
u(x 1 ,x 2 ) 



(14b) 


It should be noted that a knowledge of the function w(XpX 2 ) uniquely determines the 
functions u(x^) and u(xp,x 2 ). Substitution of equation (13) into equation (14b) results in 
the relation 



so that the function u(x^) and, hence, the function u(XpX 2 ) may be directly solved for 
once w(XpX 2 ) is known. 

Picking from a population with density function f(Xp x 2 )/w(Xp x 2 ) is equivalent to 
first choosing x^ from a population with density function f(Xj)/u(x^) and then choosing 
x 2 from a population with density function f(x 2 /x p)/u(x^, x 2 ) where the value of x^ for 
the second population has been set at the value chosen from the first population. Equa- 
tion (8) can now be written in the more revealing form, namely, 


5 



7 u(x 1 )u(x 1 ,x 2 )g(x) f(x) 


// u V 1 )»V 1 ,x 2 )(^A 1 ,, 2 )-^gL, I , i 

J- Ci 


u(x 1 )u(x 1 ,x 2 ) 


dXr 


Jj ^(x 1 )u(x 1 ,x 2 )E[g/x 1 ,x 2 l - u ( x 1 )Efg/x 1 1^ 2 ~' ( ~ f( ^l ( , X2) 


u(Xi)E[gA 


x 2 ) 


dxj dx 2 


(8a) 


In equation (8a), the variance has been broken down into components where each of the 
integrals on the right side can be interpreted in the following way: (1) the first integral 
is due to the variation in g when x-^ and x 2 are both fixed and the other variables are 
allowed to vary; (2) the second is due to the variation in E[g/xpX 2 ] when x^ alone has 
been fixed and x 2 is allowed to vary; and (3) the third integral is due to the variation in 
E[g/xj] when x^ is allowed to vary. Note that it would be an easy task to generalize 
equation (8a) to any number of sets of pertinent variables. 


FORMULATION OF THE METHOD OF IMPORTANCE SAMPLING 
FOR DISCRETE VARIABLES 

For the sake of mathematical simplicity, assume that each of the sets of random 

variables Xp . . . , x n are discrete and let us set up the foregoing problem by using a 

notation in keeping with this assumption. The results obtained can readily be modified to 

apply to the case when the variables are continuous, and the simpler discrete model will 

be used to establish the relations that arise when considering the method of splitting. 

First the notation is defined. Let the possible value sets of x^ be put into one-to- 

one correspondence with the index i=l, 2, 3, . . . ; let the possible value sets of x 2 be 

put into one-to-one correspondence with the index j=l, 2, 3, . . ., etc. Whether the total 

number of possible value sets for any variable is finite or infinite does not change the 

th 

problem. Set p. equal to the probability that x. takes on its i Ln value set and set 

^ ^ th 

p,. equal to the conditional probability that x 9 takes on its j value set given that x, 

U th ^ 1 

has its i value set. We have 
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( 1 ') 


E [g] = £ p i p ij p ijk- • • g ijk. . . = Z p i E i [g] 

i, j,k, . . . i 

E p ij p ijk p ijkT ' ' g ijkl. . . ~E p ij^ij[^^ ^ ) 

. . j 


E ij [g] = E p ijk p ijki- • • gijki. . . 

k,l,. ■ • 


( 4 ') 


’ 1] k 5 “.. 


E p ijk p ijkl. . . ( g ijkl. . . E ijf g O 


( 5 ’) 


(Primed and unprimed equations of the same number are analogous. ) 

The scheme by which the average of g is obtained is illustrated in sketch (a) where 



Pn 

u n 

fi2 Jhi, 
u 12 

P 21 n 21 | 
U 21 

P 22 n 22 
iJ 22 


911 

912 

9 2 1 

922 


the number of possible x^ and X2 (for each x^ value set) value sets are shown as two 

in the sketch. A sample size N is first pulled from the x., -population where the proba- 

til ^ 

bility of getting the i n value set of x., is p- /u.. The number in the sample possessing 

til ^ ^ ^ 

the i value set of x^ is designated as n. and the expectation value of n^ is 

N(p./u.). Then this sample is further subdivided by picking from that x 9 -population 

corresponding to the value set of x 1 where the probability of getting the j value set 

f Vi 1 

of x 9 given that the i in value set of x 1 has been chosen is p../u... (The u. and u. . 

^ 1 13 1J 1 1J 

are weight functions that play the same role as u(x^) and u^^Xg) did in the case where 

the random variables were considered as continuous. ) The number of members in the 

sample possessing the i *"* 1 value set of x^ and the j * 11 value set of X2 is denoted as 

m ^ We note that 
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(14a') 


lllllllllllllllll■llllllll 



i 


= i (i4b’) 

u ij 

3 


Let (g. .) 
value and X 2 
There are n„ 


be the a th measurement of g where x. has been fixed at its i 

.i-L 

at its j value regardless of the values of the other sets of variables, 
such measurements. Consider the random variable 


Z = 


n. 


i] 


£ £ u i u ii (s: ii> 0 

i,j a=l 


N 


(16) 


Equation (16) can be written 


Z = 


£ z « 


where 


u.u. . 


n ij 


z ii - E ( ^ii) 

13 N 13 


a 


0=1 


u.u. . 


= - — - 3 n- -(g- •) 
N 13^3 ; 


Taking the expectation value of Z. ^ gives 


(17) 


u.u. . 

E[ z ii] =- LJU E 
13 N 


- — 

u.u. . 

- 

1 - -|“| 

n ij (g ii>_ 

- 1 1] E 

N 

V 

(g. .)/n. . 

L 13 13 J 


8 


where E 


(g. . )/n . . 
B i]' iJ 


is the expectation of the average measurement of g.^ when the num- 


ber of measurements is n^. Because this expectation value is independent of the num- 
ber of measurements, the equation becomes 


E[Z,.] = E [g]N P J. -U = PP] E [gj (18) 

13 N 13 u t u t - 1 13 13 

so that 

E[Z] = Z E t Z ij! = Z p i p ij E ij[sJ = E [g] (19) 

ij i, j 


Hence Z is an unbiased estimator of E[gl. 

9 

We wish to find the variance of Z. The following relation holds: 


4 = E t z2 ] ' e2 [ Z ] = Z E f Z ij Z i'j-] - E 2 f g ] 


i,J, i’J’ 


= E E [ z fJ + E E [ z ,j z ii'i + E E t z ii z ri- 

i, 3 i, 3, 3* i,3,i',3’ 


- E 2 [g] 


(20) 


(3^3’) 


(i^i') 


Considering the terms of the first summation on the right side of equation (20) gives 




2 2 


2 2 

_ 


E 

cm ; 
lN 

. U l^l E 

2 / — \ 2 
n. .(g, ■) 

= Til E 

2 rp 

n. .E 

(g.,) 2 /n. • 


L 13 J 

N 2 

L 13 13 J 

N 2 

13 

_ 13 l 3j 


2 2 

u A>e 

N 2 


n £ ^ + E fjg] 


13 L 


2 2 
u. u. • 

1 13 


N 


2 2 
u, u. • 


< E[n. j ] + J-llEf.[g]E 


g 


13 


N 


13 L 


n 2 . 

13 


The relation 
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E["ij] = - — 

1 u i u ii 


has already been employed in obtaining equation (18). To find an expression of E 


2 

n . . 

i] 


we employ the fact that, for a specified value of m, the quantity n^ is distributed as a 
binomial variable with probability p.j/u.^ of "success" where the number of trials is 
equal to m. Hence 


E 

r 2 

n. . 

= E 

E 

[2/1 
n. ./n. 


_ i] 



L u i \ 


= ^ii E 

r 2 " 

n. 



2 

U ij 

l 

u.. 1 
i] 

V “ J 


2 

2 
u. . 
i] 


Pi / Pi' 


N 2 ^ +Nli[l-- 1 
2 u. \ u. 

“i i \ V 


P i P ii I P il 
+ N — 1 - 1] 


u. u. . 
1 IJ 


u. . 
1] 


o P? Pfi Pi Pii / Pi Pii 

= N Z -L -il + N -i fl - -i -H 

2 2 u. u. . \ u. u. . 

u i u ij 1 U V 1 


Thus, the terms of the first summation are given by 


*3 


2 2 r 2 r .1 „ z l '•"l * 11 1 1 1 ii i ^ ^ r ^ 

p. p. . E . .] g + — p.p..u.u. a + — J 1 1 U-U..E.. g 

1 N 1 y 1 H Sii N u u ' -- - I 1 H J 

J i ij 


Pi Pii / Pi Pii \ 2 2 e 2, 

ij ij 


u. u.. 

1 H, 


( 21 ) 


Considering the terms of the second summation gives, for j^j', 


E[Z Z ] = TVt)’ E lj [g|E ij ,lg]E[. W | 
1ST 


where we have used the fact that (g^) and (g.^,) involve measurements made on two 


10 



different groups of the sample and, hence, are independent variables. We have 


E[n..n..,] = E n..E[n. .,/n. .] 

1 i] j L i] L ij ir 


Now with n,, fixed, n.., is a binomial variable, where the number of trials is (n. - n. .) 
i] 5 i] i 

and the probability per trial of taking on the j' value set of x 9 is (p..,/u..,)/ 

U Ij 1 J 

[1 - (p../u..)]. The denominator 1 - (p../u^) is needed as a normalization factor modi- 
fying the probability p..,/u. because, with n.. fixed, none of the (n. - n..) trials can 
. . ^3 ^3 fjn ^ 

give rise to a member having the j value set of Xg. Thus 


E rv /n ij] 


Pii- 

<"t - f- 

ij’ 

i.!a 

u ii 


and 


E Kj n ij'! 


U ij\ 


111 

U ij 


n.n. . - n. . 
i i] ij 


Continuing, 


E [ n i} n ij' ] 


U ij T 


E 


r 

p< 


p. . 

1-111 

U ij 


11 n 2 „ / n 2 . 121 h - 


u. . 
13 


i + 



p yj p ii 

U ij’ U ij 


2 

n i " n i 


_ P U’ P ij P i 

V U i3 uf 


N 2 fl - - 
N 


11 



Hence, for j/j’, 


E I Z ljV = p l p iJ p iJ ,E ijW E ij , W ^ - 5) 


(22) 


Considering the terms of the third summation of equation (20), we have, for i/i’. 


E [ Z ij Z i T ) =— f 1 ' 1 ’ E Ij WE 1 . 1 .te]E[» ij n 1 . j .! 


where 


E [ n i) n i'i , l = E [ n ij E [ n l'j' /n iji 


- E 


Pi* Pop 
(N - n..) -1 _LL 




n ij 


u i’ U i’j’ 


p. p.. 
1 _ 1! 
u. u.. 

1 13 


P_v Prr 

U i’ U i’3- E 


1 - 


p i Pii 

u i u ij 


Nn. . - n - 
i] 13 


Pi, 

Pi* j' p i 


u. f 

U-, 

u. 

u, ■ 

x T 

1 3 

1 

13 


N, 


Hence, for i/i', 


E[Z i jZ iljI ]=P i P i .P i , Pi , if E i j[gJ E 1 . j .fe](l-^ 


(23) 


Substituting equations (21) to (23) into equation (20) and using the relations 


V 8 1 = % + E i/ g ' 


O') 


and 


w. . = u.u. . 
ij 1 13 


( 13 ’) 


12 



yield 



X 

i,3 



Pi Pi 


1 ? (’ - ? ? )Ww 


u. u. . 

1 l], 


- z 

i, j, j 
(j/i') 


Pi p ij p ij’V g ] E ij ,[Bl 


Z PiPi'Pi3 P i’j- E ij[ g ] E i’j’[ g ] 
(i/i’) 


'IwL + 2 iVij u i u ij E «[ g ] - 



= Z w ij E i/ g2 i p i p ij - e2 ( s| (8,) 

i, j 


Equation (8'), which is analogous to equation (8), can be written in a form analogous to 
equation (8a) as 


No, 


Z 


2 2 
w. .a 
i] g. 


P i P ij 


i] w. . 
i] 


L 


EJ 


w. .E. .[g| 
1} iJ l6J 


u.E.fg 

l i Lto 


2 p i p ij 


w 


i] 


Yj^W - E[gJ ) 2 ~ 


The analogous equations to equations (10), (11), and (15) are 


w 


ij 



(8a’) 


(10’) 


13 




MODIFICATION TO INCLUDE THE METHOD OF SPLITTING 


We now modify the scheme by which the average of g is obtained in accordance 
with sketch (b). 



P 11 


1 

W U _ 


U 11 


S 11 


911 

_^12 

p n 

1 

u 12 


u 12 


s 12 


912 

£21 

n 21 

1 

U 21 , 


U 21 


S 21 


921 

p 22 

n 22 

1 

u 22 , 


u 22 


s 22 


9 22 


(b) 


Again the possible x^ and X2 (for each x^ value set) value sets are shown in sketch (b) 
as two in number. 

The first step in the sampling procedure remains the same; a sample of size N is 

iL 

first picked from the x 1 -population where the probability of getting the i 111 value set of 
x 1 is p./u^. The number in the sample possessing the i Ln value set of x^ is desig- 
nated as n^ and, for a sample of size N, the expectation value of m equals (p./m)N. 
This number n. is then multiplied by the splitting factor l/s- to yield the number 
ia - (l/s-)n i that now possesses the i Ln value set of x^. (It may be noted at this point 
that no longer does the sample necessarily consist of N members because ^ ia ^ N 

i 

in general. ) The sample is further subdivided on the basis of the variable x 9 and n.. 

th ^ 

represents the number of members in the sample with the i value set of x^ and the 
j tJl value set of Xg where the probability that a member possessing the i** 1 value set 


14 



I 


th 

of x, also has the i value set of x 9 is p,,/u. .. Again, the splitting factor l/s,. is 

1 “ *1 th l J 

introduced and the number of members possessing the i value set of x. and the 
j value set of x„ is changed to v.. = n../s... Finally, measurements of g.. 
are made where the a measurement is denoted as j)^- 

It is observed that two steps are included in each stage: (1) importance sampling 
where the weight factor l/u. or 1/u.j alters the selection probabilities and (2) splitting 
where the splitting factor l/s^ or l/s^ alters the numbers selected. Such a stage will 
be designated as a ''composite” stage. In general, whenever a splitting step is present, 
the total number of members in the sample changes. 

The random variable that is of interest now is 


"ij 


> > u.s.u. .s. (g- .) 

/Li sLi ill] ij v& lj' 


a 


g = i, 1 «=1 


N 


(24) 


Proceeding in the same way as before, we write 


3 



3 ij’ 3 



u.s.u. .s. . v. .(g. .) 

l l i] i] i] v& ij' 
N 


(25) 


E 13 


i]' 


u.s.u. .s. • 

1 LiiJi E 
N 


v. . 
i] 


Z <gn> 

a 1 


i yot 


u.s.u. .s. . 

= 1 1 l l U 
N 


EijlgM^ij) =P i P ii E ii [g] 


1J 


i] 1J U 


(26) 


E iji = Z E i/i j i=Z p i p ii E ii [B| - Elg| 

i,i M 


Hence, is an unbiased estimator of E[g], Continuing, we have 


(27) 


15 



( 28 ) 


»} =E[J> 2 |- E 2 [J>| 

= Z e ^]+Z e! ^> + Z - e2[sj 

i,j 

(i/r) (i/i’) 


As before, each of the summation terms of the right side of equation (28) will be 
considered in turn. For the first summation terms we get 
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This yields 
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For the second summation j/j' and 
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and we use the method of Lagrangian multipliers. Equations (34a) and (34c) merely re- 
quire that the altered probability distributions be true probability distributions whereas 
equations (34b) and (34d) require that, on the average, the number of members in the 
sample after the splitting is changed by the factor for the x^-variable stage and 

mg for the Xg-variable stage. (From this point of view, the factors m^ and mg can 
be regarded as "magnification" factors. ) 

Multiplying each of the equations of constraints by so-called Lagrangian multipliers 
and adding these terms to the expression for Ncto defines the quantity 


L = NoJ, + A 




+ A 


III 


PiPij 

W ij 


- m , m 


12 


L 1 


where Aj, Aj p A^(i=l, 2, . . . ), A pi denote the Lagrangian multipliers. This expres- 
sion is minimized with respect to each of the four variables u p w., w.j, and w^ by 
setting the appropriate partial derivatives to zero: 
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Solving equations (34) and (35) for the optimum values of m, w^, u^, and gives 


A E iM 

Sj E fei 
JL = E ijW 


_J_ = 1 

2 g 

y 

w.. E 
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(36a) 


(36b) 


(36c) 


If m, m., and w^.. take on the values given by equations (36), then it turns out that 
NcJ? is independent of the form of w^ so that the simplest expression for l/s^ in ac- 
cordance with the constraint given by equation (34b) is 



which yields 

! mjEjg] 

Wi " Efg] 


(36d) 


It can be seen from equations (36a) and (36b) that the optimum values m and u^ of 
the weight factors correspond to ’’adjoint biasing" (refs. 1, 10, and 11) wherein the 
biasing as given by the reciprocal of the u's is proportional to the expected contribution 
of the member to the answer E[g]. 

The minimum value of Nero is found by substituting equations (36) into equation 
(32a) " 


XT 2 

NC 2 


v. 


m l m 2 


(37) 
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The previous equations can easily be converted to the case wherein and x 2 are 
continuous variables (or sets of continuous variables). The following relations hold 


u(x 1 )s(x 1 )u(x 1 , x 2 )s(x 1? x 2 ) ^g(x x , x 2 )^ 

3 = — 


a. 


N 


where the summation is taken over all measurements. 


E[J] = f J E[g/x 1 ,x 2 ]f(x 1 ,x 2 )dx 1 dx 2 - E[g] 

No) = Jj w(x 1) x 2 )(a|/x r x 2 )f(x 1 ,x 2 )dx 1 dx 2 

+ jj w(x 1 )u(x 1 ,x 2 )E 2 [g/x 1 ,x 2 |f(x 1 ,x 2 )dx 1 dx 2 
- J w(x 1 )E 2 [g/x 1 ]f(x 1 )dx 1 
+ J u(x 1 )E 2 [g/x 1 |f(x 1 )dx 1 - E 2 [g] 


w(Xj) s u(x 1 )s(x 1 ) 

w(x 1 ,x 2 ) H u(x 1 )s(x 1 )u(x 1 ,x 2 )s(x 1 ,x 2 ) = w(x 1 )u(x 1 ,x 2 )s(x p x 2 ) 
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(34c’) 
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CARRYING OUT OF THE SPLITTING PROCESS 


The question now is how best to carry out the process of splitting. According to 


sketch (b), in the splitting part of the jx^ s ^ a S e > the number j n *J> rnern t >ers is mul- 


f 1/s i 1 K 1 

tiplied by the factor j^y s J to yield the number | • This is a simplified concept 


of 
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what can be done because ^iy s j is not, in general, an integer. Since, in a literal 

sense, we cannot work with fractional members, a method must be devised to give the 
necessary flexibility. 

In some cases, the number / 1 of members corresponding to a given value of 


n. 


n * j is left to chance where the probability is so chosen that the required relations 


Efi/./m] = — n. 


(38a) 


s i 


E[v-./n. ■ ] = — n. . 
L i] i] ' e i] 

S ij 


(38b) 


hold. For example, suppose we write 


s i 


= L + n 


(39a) 


s, j 


= I.. + r.. 
13 i] 


i] 


(39b) 
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r j 1. Then, if each of the 


where ( j 1 > is a non-negative integer and 0 

l ijJ 

fL ] fl-r. ] r i i +1 l f r i] 

members gives rise to < j 1 > members < , 1 f of the time and < j - > members < 

l ijJ l ijJ l ij J l ijJ 

fn. "j fl.+r. ") 

of the time, on the average, each of the J 1 l members generates -j j 1 * f members, 

t ijJ l ij ij J 

thereby satisfying equations (39). This process is an example of the "Russian roulette" 
method wherein particles are created or annihilated by chance. If this process is used 
as described to obtain the sampling as a function of the variables x^ and Xg, then the 
expression for the variance of £ as given by equation (32) no longer holds. The reason 

f v i ] 

is that additional uncertainty has been introduced by the fact that the variance of J V 

l ijJ 
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corresponding to a given value of i * i is no longer zero. For this technique, the num- 

l ijj 

fn. •) fL+n 

ber of times 3T that the | 1 l members give rise to -s j 1 + J members is a binomial 

l ijj l ij J 

variable where ft) is the probability of "success" and ft.) is the number of trials. 


*3. 
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it is easy to show that 
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1 ' 1 
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and 
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2 +n tj r ii ( ‘- r ii ) 


(40b) 


With the aid of these relations, the expression for the variance may be obtained in 
exactly the same way as carried through in the preceding section. The random variable 
£ is given by equation (24) and each term of equation (28) must be reevaluated for the 
situation under consideration. A point that might cause some difficulty is the evaluation 
of E[ia ] for j^j’. This is accomplished as follows: 

E ( Vir ' = ffffff V W hi- V^i d h dn u % • d hj d V 


where f(n., ia, n.^, m.,, y m.,) is the appropriate multivariant probability density func- 
tion. This equation can be written as 
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E [Vij'l = iff f lVrVV l *i dl 'i dn ii %■ 
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/ f ; V ! 'ir f ' : i;- V iy /n v "l- V d V 


But because the random variables ia • and only depend on n.j and n.^,, respec- 

tively, we get 


E[ Vir 1 = ffff St "> V V n ij’ )d "i d "i dn ii dll ii' 
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For the splitting technique involving Russian roulette as described, the increase 
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in E 
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over that value given by equation (29) is found to be 


A E 


> 2 . 


I PiPj^iVid - r i) E f 3 [e] 


+ 5 PiPij w ii s ii r ijd - V E ^ gl 


The increase A (e ,]) in E J.., | over that value given by equation (30) is 


found to be 
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Equation (31) is found to remain unchanged. The increases result in an increase 
A^N in NcrJ> over that value given by equation (32): 

A ( Na !) = £ p i w i s i r i (i - r i )E x 2 [g] + £ p i p ij w ij s ij r ij (i - r ij )E i 2 jt g i (4D 

i ij 


To get an approximation of the magnitude of this increase, we make the tentative 
assumption that both m and r.^ are uniformly distributed over the interval (0, 1) so 


that 


E[r( 1 - r)] 


* / n 


r(l - r)dr = 


(42) 


and equation (41) becomes 


A ( Na J) - eZ PiWiS ^ [gl + ^Z PiPijWijSijE ^ [sl 




(43) 


In general, the increase in the variance of £ as given by this equation is not insig- 
nificant and, in some instances, could practically nullify the reduction obtained by split- 
ting. Hence, the method of Russian roulette should be used with caution. 

A technique v hereby equations (38) are satisfied without the introduction of additional 

f n. i 

uncertainty in the final result is that in which, for each of the original members, 

I / +if members are generated, the last member having a weight 1 r * > that of the first 

l ij J l ijj 

fL 'j fr. "I fl.+l "i 

J members. Hence, in this method, provided ^r..j ^ 0> we actually follow 

members but weight the last one differently. Here our relevant random variable has 
changed from that given by equation (24) to 


bi 


(44) 
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where 
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(45) 
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and we have subdivided the paths as shown in sketch (c). Hence 
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(46b) 
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Note that 


w. . 


KIV ij. = N li E ijfe]N hhilVij + r i I ij + ■'ij'i - r i r ij) = "u E iifel a Z*!z = >VV' ijl> 1 


Pipi 

i u ij s i~u 


and thus 


E[Y] = E[g] 


as desired. 

This choice of random variable introduces no new uncertainties into the calcula- 
tions and the fact that we are dealing with more members than previously (because 

) results in a reduction in the variance. In appendix A, 





it is shown that 



2 

where No^ is given by equations (32). The disadvantage connected with this technique 

is that only the fractions r. and r.. of the measurements made on the last members of 

1 i] 

each stage are being used. In this sense, the method is not as efficient as it might be. 

A more simple way of proceeding is to insist that l/s^ and l/s^ are both integral 
for every i and j. In other words, even though the optimum values of these quantities 
are nonintegral, we always take, as the value to use for l/s^ or the smallest 

integer larger than or equal to the optimum value. In this way, we avoid both of the dis- 
advantages associated with the two methods previously described - namely, (1) no new 
uncertainties arise and (2) each of the (g^)^ measurements corresponding to particular 
values of i and j carries the same weight. However, if this method is used, then, as 
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also occurs in the previous method, the number of particles in a splitting step never de- 
creases. 

Kahn (ref. 3) suggests two separate treatments in a splitting step and which of the 
two treatments an individual member receives depends on its values of the measured 
x-variables. Sample members in a splitting step are classed as belonging to either a 
type I or II region. In type I regions, the optimum value of l/s^ or l/s^ is less than 
or equal to, say, 0. 5 so it is desired to decrease the number sampled in these regions; 
in type II regions, the optimum value of l/s. or l/s^ is larger than 0. 5 so the num- 
ber sampled in these regions should not decrease. Kahn uses Russian roulette on type I 


region members whereby a member having an optimum value of 


0. 5 is given the chance of 



of going on and the chance of 



less than 


of being killed, 


thus satisfying equations (38). For type II region members, "integral" splitting is 


used where the 


value of 



taken on is the reciprocal of the integer closest to the optimum 


After going through two stages of splitting in each of which a mem- 


ber can be classified into one of the two regions I and II, there are four classes of mem- 
bers, that is, I- 1, I-II, H-I, and II-II. Equations (32) and (41) apply in this case to 
yield the proper expression for the variance. 


EXTENSION OF ANALYSIS TO ANY NUMBER OF STAGES 


Equation (32b) gives the expression for variance where two stages have been em- 
ployed, each stage being a composite stage consisting of an importance sampling step 
followed by a splitting step. This equation may be easily generalized to apply to any 
number of stages. For example, if there were three composite stages, then the appro- 
priate relation becomes 


No 2 = 
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(48a) 
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where 


w ijk ~ w ij u ijk s ijk 


W ij = W i U ij S ij 


w i = u i s i 


(48b) 


and where the symbol p has been introduced to denote the unconditional probabilities - 
namely, 


/&i = Pi = probability that takes on its i** 1 value set 

th 

p = P^P^j = probability that takes on its 1 value set 
th 

and X 2 its j value set, simultaneously 

^ijk = ^ij p ijk = p i p ij p ijk = P robabilit y that x i takes on its 

i^ 1 value set, Xg its j tfl value set, and Xg its k^ 
value set, simultaneously 


(48c) 


In general, for any given number k of successive composite stages, the equation 
for variance can be written 
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or 



J.L 

where is the probability that takes on its iy 1 value set, p^ . is the probability 
1 12 

that both takes on its i^ 1 value set, and Xg its i^ 1 value set, etc. 

Again, as noted just after equations (33), these expressions for the variance also 
hold for the cases where the stages are pure importance sampling stages (s=l) or pure 
splitting stages (u=l). (We shall designate a composite stage by the symbol "4'", a pure 
importance stage by "I”, and a pure splitting stage by "S".) Also, as demonstrated 
previously, these expressions may be readily modified to apply to continuous variables. 
Finally, these expressions must be adjusted in accordance with the relations of the pre- 
ceding section if applicable. 


OPTIMIZATION OF THE WEIGHT FACTORS FOR NON-NEGATIVE g 

If, for the moment, we ignore how the actual splitting process at any given stage is 
to be effected, expressions for optimum weight factors may be worked out. For example, 
in the case of two pure importance stages (I-I), the weight factors u^ and u^ satisfy 
equations (14a') and (14b'), respectively, and the optimum choice of these weight factors 
is given by equations (10’) and (15') where equation (13') applies. For the case of two 
composite stages, 'f r -4 r where g is non-negative, the weight factors satisfy the condi- 
tions of equations (34) and their optimum choice is governed by equations (36). 

For the general situation, it should first be noted that the number of stages in the 
sampling procedure does not necessarily change the minimum variance attainable. As an 
example, let us again restrict ourselves to non-negative g and consider the three stage 
sampling I-S-4>, that is, where the first stage is a pure importance stage, the second 
pure splitting, and the last a composite stage. In this case, Spl and u„=l so that 
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W ijk = u i s ij u ijk s ijk J 

and equation (48a) can be written 

N “ 2 = X ' /L ^i]k W ii U ijk E l 2 jk[®J 

i,j,k i,j,k 

- 2^ij W ii E i 2 jf g] + Xl ^ij U i E ii [g] ' E2[g] 
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The constraints to be satisfied by the weight factors are 
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m 2 


(52a) 


(52b) 



i, i,k 


(52c) 


(52d) 
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Using the method of Langrangian multipliers to determine the minimum of Na sub- 
ject to the constraints of equations (52), we find that the reciprocals of the optimum 
weight factors satisfy 
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where is defined by the equation 
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It turns out that the minimum value of Na is independent of the form of 
equal to 


NT 2 

Na = • 


E 2 

k 1 


L & ijk 


m 2 m 3 


+ E>j] - E Z [g] 


where the sum of the last two terms on the right is non-negative as shown by 

E 2 ^] - E 2 [g] = (E[e,] + E[g]) • (E[ fi ] - E[gj) 

= (e[ c j] + E[g]) • (e[£j - EjIgJ) 

coupled with the relation 


4 - Ef[g] = E i 


( E ij[g] " E i[g]) 


Another way of demonstrating that ^ E^[g] is 
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However, if the splitting stage in I-S-'k is eliminated by changing to an 
sampling where ( '^) 2 3 denotes a composite stage for the two sets of variables Xg and 
Xg simultaneously, then, equation (32) applies where j — (j,k) so that 


No 2 - ^ i“iOk) w i (jk)°g 1(Jk) + X Pl (jB u i u i(ik) E f(jk)W - E 2 fe] 


i, (j,k) 
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The constraints to be satisfied by the weight factors are 
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The reciprocals of the optimum weight factors for non-negative g are given by 
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and the minimum variance for non- negative g is 


/X 

No 2 = 
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g i(jk)_ 


(59) 
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23 


Hence, if nigg of equation (59) is equal to nigm^ of equation (55), then the minimum 
variance of an sampling is less than or equal to the minimum variance of an 

I-S-'k sampling. (Of course, in order to attain the minimum variance in both cases, 
the answer E[g] among other things must be known before sampling begins. ) Note that 
the optimum biasing of the I stage in the sampling is adjoint biasing whereas, 

in the I-S-'k sampling, the optimum biasing of the I stage does not correspond to ad- 
joint biasing. 

Another point may be discussed in connection with our problem. If we consider a 
('T')l 23 sampling, then the minimum variance expression is again given by equation (59) 
with n^g going over to m^g and the optimum biasing of the importance step of the 
composite W^g stage corresponds to adjoint biasing, namely, 


U (ijk) 


li)k) 

E[g] 


[g] 


Again the minimum variance attainable is not affected adversely. 

It must be mentioned here that a three- variable sampling yields an optimum 

minimum variance that is, in general, less than a single-variable '!> sampling (contain 
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ing the same number of elements in the sample) wherein only the first variable x^ is 
measured. For a sampling and non-negative g, the minimum variance is 


/s 

NCX 2 
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(60) 


whereas equation (59) with m.^ replacing irigg holds for the (T) ^23 sampling. Inas- 
much as 
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cr 

> E 
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1 gi J 
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as demonstrated in appendix B, our contention is readily established. 

For non-negative g, it turns out that with any given number k of sampling stages, 
the minimum variance with the least number of members in the sample is obtained if the 
last stage is composite and all stages except the last are adjoint-biased pure importance 
sampling stages. The importance step of the composite stage should be adjoint biased 
and the splitting step should be biased in accordance with 

m a 

■ 

Tw ■ 

The minimum variance for this situation is 
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(63) 


It is instructive to investigate the general problem of the optimum biasing of a stage 
for other cases than that just presented. It is found that, in order to set up the optimum 
biasing of a stage, one must know the nature of the stages following the stage in question. 
The optimum biasing results for non-negative g are depicted in the accompanying 
table I. The following notation is used in the table: 
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5 2 - e 2 - E 2 [g]; 5 2 , e 2 - E 2 [g]; 6 2 , e 2 - E 2 [g]; 


y 2 (D) S E 2 [D] - E 2 [g] 
y 2 0^E 2 0-E 2 [g] 
^(□)-E 2 [n]-E 2 [g] 
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where the symbol [J stands for any quantity. For example, by equations (67), 

nj< £ ijki> = E fj[ £ iiki] - E i 2 jle] 

= E . 2 [ T ij] - E f[eJ 


The symbol m. is the magnification factor for the i^ h stage. Of course, mu = 1 if the 

iu 1 

i u stage is a pure importance sampling stage. 

Equation (Bla) shows that 


and equation (Bib) shows that 
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so that 
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( 70 ) 


T i = E i 


‘ 2 


2 

T ijkj 

= E i 

T ijkl . . . 


As can be seen from table I, the optimum importance sampling step of a composite 
stage is always adjoint biased. Also, in general, the optimum choice of a weight factor 
for a particular stage depends on the nature of the stages following the stage under con- 
sideration - in particular, on whether the stages terminate before the first pure splitting 
or composite stage occurs. (Note that the results are the same for the pure splitting 
stage and the splitting step of the composite stage. ) The blank spaces in the tabulation 
for the S or steps indicate that the minimum variance is independent of the biasing 
employed with these stages and, hence, for simplicity, uniform biasing is to be em- 
ployed. 

As an illustration of the use of table I, let us consider a four-variable sampling of 
non-negative g which is to proceed as 'k-I-S-I. The following conditions apply in gen- 
eral to the weight factors: 



According to table I, the optimum weight factors are given by 
1 st stage: 'J> followed by I-S . . . 
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Importance step 1^: 


i = w 

Gj E M 


Splitting step Sy 


1 


w i 


m in (e ij> 

e ivy 


2 nc * stage: I followed by S . . . 


e. . 
i] 


u « 


E i[ e ijl 


r H 

3 stage: S followed by I 

= m i m 3njk< T i3ki> 
w ijk E hjk^ T ijkl^] 

4 th stage: I (last stage) 


1 . T ijkl 

i;kl E ‘iK [T liM 1 


(It may be noted that eq. (72e) arises from an "extrapolation" of table I. ) 

2 

these values into the expression for Nct yields the result 

T 2 . + 

m^N m^m^N 


where each of the terms on the right side is non-negative. 


(72a) 


(72b) 


(72c) 


(72d) 


(72e) 

Substituting 


( 73 ) 
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OPTIMIZATION OF THE SIZE OF THE SAMPLE AT EACH STAGE 


The magnification factors and mg may themselves be optimized. Suppose, on 
the average, the total cost T of conducting the Monte Carlo analysis for the example 
I'-I-S-I just considered is given by an expression of the form 

T = c^N + CjgtnijN) + Cgg(nijN) + Cj^mjmgN) + c g (m 1 mgN) 


where 

c^ = average cost per sample member processed through composite stage 1 

C I2 = average cost per sample member processed through importance stage 2 

C S3 = avera § e cos t per sample member processed through splitting stage 3 

C I4 = avera S e cost per sample member processed through importance stage 4 

c = average cost per sample member of measuring g 

o 

The above equation can be written as 

T = a^N + a^m^N + a^m^rngN (74) 

where 

a 2 - c 
a l - c *i 

a 2 = C I2 + C S3 
a 3 = C I4 + c g 

In practice, c^ and Cgg are not constants independent of m^ and mg, respectively, 

but in this analysis we shall assume that such dependence can be neglected without in- 

/\ 

2 

troducing appreciable error. The problem as now set up is the minimization of a 
given by equation (73) with respect to the variables N, m^N, and m^nigN subject to the 
constraint that T of equation (74) is fixed. Carrying out this minimization process re- 
sults in the theoretically optimum values 
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N = 0 




TE 


(miN) = — r — V 


[ y i (e ij>] /a 2 


1# a 2 E Mij)] +a 3 E hk< T ijkl)] 
TE [ y iik^ T iikl^l / a 3 

(m m N) = k f - 

a 2 E [^i^iP] + a 3 E [VW] 

and the minimum variance for given T is 

+a 3 E hk< T iiH)] 






(75) 


(76) 


where the corrections to be made dependent on how the splitting processes are carried 
out are not included in equation (76). Of course, since N must be a positive integer, 

the theoretically optimum result N - 0 is not allowed and, in actuality, N should equal 

A 

the smallest positive integer, namely, 1. The reason that N turns out to be zero is 
that the true optimum values of the weight factors are being employed. Appendix C 
shows that any deviations of the weight factors from their optimum values result in an 

increase in N. 

It was remarked previously that the corrections dependent on how the splitting 
process is effected are not considered in equation (76). As an example of how such cor- 
rections may be included, we turn again to the three-stage sampling I-S-’k considered 
previously. The minimum value of the variance is given by equation (55) as 


2 


a 


y 2 ^) 

+ 

N 


E' 





(55a) 


where the notation of equations (65) and (67) had been incorporated in equation (55a). 

Here we also wish to obtain optimum values of N, n^N, and n^m^N subject to the con- 
straint that the quantity 


T = bj N + bgmgN + bgm 2 m 3 N (77) 

It is realized that equation (51), from which we obtain equation (55a), does not include 
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the additional variance terms (such as those given in eq. (41)) for the S and 'I' stages. 
Consequently, equation (55a) cannot represent a true minimum. Let us first consider 
the S stage. It is assumed that the added terms, that is, the second summation of 
equation (41), constitute only a small perturbation. If this is so, then equations (53a), 


2 

(53b), and (53c) and the result that a is independent of the form of w^ and, hence, 
on the form of s^ can be considered as approximately correct. For simplicity, we can 
take s^ to be a constant independent of i and j where, by equations (52a) and (52b), 
this constant must be l/nig. If Kahn's procedure is followed, then mg is either a 
positive integer (corresponding to a type I region) or mg is a fraction less than 1 (cor- 
responding to a type II region). In the first instance, mg is a positive integer and all 
r„ are zero so that the added terms of equation (41) vanish and equation (55a) is exact. 

2 2 
Now, since ct is independent of (nigN), if the process of minimization of ct were 

mechanically carried out, the optimum value of (mgN) would be zero. But, because mg 
is restricted to the positive integers, the lowest value for mg that can be chosen is 
unity and the splitting stage S is reduced to a unit stage U where s^- = 1. In the sec- 
ond instance, mg is a fraction less than 1 and equation (51) must be modified to include 
the appropriate terms of equation (41). Minimization of the variance with such terms 
included show that the optimum value of mg for this case is unity. Thus, it has been 
proved that, to obtain the minimum variance, the splitting stage becomes a unit stage 
whereby the three-stage sampling process goes over to I-U-'k where equation (55a) is 
now 


2 


CT 


y 2 ^) 

N 




m 3 N 


(78) 


Equation (77) becomes 


T = bT + b« N 


+ bglUgN 


(79) 


(It should be remarked that this result holds in general; that is, when the variance 
of the sampling process does not depend on the form of the splitting factor for a particu- 
lar S stage or for the splitting step of a composite stage (blank spaces in table I), in 
order to obtain minimum variance, the splitting factor goes to unity. ) Minimization of 
2 

ct of equation (78) with respect to N and m^N subject to T of equation (79) being 
fixed yields 
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and 




(80) 


( 81 ) 


It must be remarked that the additional variance due to the splitting step of the <k stage 
which is equal to 

X > s ijk w ijk s ijk r iik< 1 - r iik> E ?jk! g ! 
i, j,k 


has not been incorporated into equation (81) and the optimum values of 
be readjusted to keep the additional variance small. 



may have to 


REMOVAL OF NON-NEGATIVE RESTRICTION 


Heretofore we have limited g to non-negative values. The hypothesis that g be 
non-negative meant that all expectation values of g are necessarily non- negative and 
permitted us to express the weight factors, which must themselves be non-negative 
quantities, in terms of these expectation values. Table I holds only for the case where g 
is non-negative. 

If g can take on negative as well as positive values, then equations such as (36a) 
and (36b) cannot, in general, be written. To demonstrate the changes in the relations, 
let us treat the same problem (a 'k-ik sampling) which led to equations (36a) and (36b) 


but now no longer regard g as being restricted to non-negative values. 

AAA A 

(34), and (35) still apply and the optimum values of u^, w^, u^, and w„ 


Equations (32), 
are given by 
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1 °T 


Ui E [ a i] 


i_ = m in (ff i[) 


w. E 


a. . 

i] 


~ Ejo'.-l 
u.. i L 11 J 
i] J 


1 2 g. . 


W ii E 


g ij 


(82a) 


(82b) 


(82c) 


(82d) 


where oo denotes the absolute value of Ej[g] and a?_ the absolute value of E-[g]; 
that is, 


a. = 

1 


E i[g] 


: 01 .. u 

’ 13 


E ij[g] 


(83) 


Note that 


a. = 

1 


EJg] 


E p ijV g| 


E p ij 

3 


(84) 


so that 


a. < EJo;..! < E.ra.., 1 — 

1 i L ij J i L ijk J 


(85a) 


Similarly, 


a.. < E.-Ta.., 1 < E..ro'.. 1 ,1 < . . 

H ij L ijk J ij l i]kl J 


(85b) 


2 

The minimum value of NuJ) is found by substituting equations (82) into equation (32a): 
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( 86 ) 



E 


g. 


i] 


E 


n< a «> 


m l m 2 


m. 


+ r 2 ( Q ' i ) 


Table II is a listing of the optimum values of the reciprocal weight factors when g 
can take on negative as well as positive values. 


A SIMPLE EXAMPLE 

We take a simple example to illustrate the concepts. We suppose the situation as 
shown in sketch (d). There are two stages of sampling, the parameters being given in 



P U = 0.90 

— EnM-o. o gii =ioo 

?\2 ' 02 

— * E 12 [g]- 100. o gi2 -50 

Pl3 ~ 0. 08 

E 13 M = 6 °g B ■' 2 

p 2 j = 0 . 20 

— - E 2 i[g] = 10 . 0^=20 

P 22 “ 0. 50 

— - E 22 tg] = 60. o g22 = 10 

P 23 ~ 0. 30 

" E a [ 3] = 0 , o g23 - 0 


Id) 


the sketch. The formula for the variance for the two-stage sampling is given as equa- 
tion (32a) to which must be added equation (41) to take into account the added variance 
introduced by the splitting processes. 

The variance resulting from straightforward sampling (U-U) for this example is 




- E 2 [g] = 3605 


(87) 


2 

We shall determine the minimum variance a for each of the following four cases: 

(1) U-'k, (2) I-I, (3) 'k-I, and (4) 'k-'k. The magnification factor of any splitting stage or 
step will be set close to unity. 
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I 


Case (1) U-* 


Equation (32a) becomes 


a . .w. .a + 

/ i] gy 


y^..u..E 2 

/ i] i] 


[g] - E^g] 


From table I, the optimum choices of the weight factors are 


i_ m2 % 

“ii E ‘ lB ’ ®u E p Bij 


Substituting these values into equation (88) results in 


E 2 <r 


— + 6 2 = 1132 + 183 


where a correction based on the splitting step is still to be added. 

For m 2 = 1, the reduction in variance over that of the straightforward sampling 
situation can be written as 


//\ 

XT / 2 2 
N CT UU <T U'J> 


m 2=! 


E r 2 - E 2 a_ - e 2 


= E (a - E a 

l g ii g ii 


+ E E [g]- E.[g] 


= 1673 + 616 = 2289 


Calculating the optimum weight factors for m 2 = 1 from equations (89) and the 
relation 


I 



1 


gives 


i 

j 

i^i 


1/S ij 

1 

1 

0 + 

2.97 

oo 


2 

40. 32 

1. 485 

0. 0368 


3 

2. 42 

0. 0594 

0. 0245 

2 

1 

0. 3126 

0. 594 

1. 90 


2 

1. 875 

0. 294 

0. 157 


3 

0 

0 



Following Kahn's procedure as described previously, we take l/s.^ and l/sg^ to be 
integers. In particular, l/s 2 ^ is taken as 2. This changes the numbers for i=2 and 
j = l from those given previously to 

= 2, — = — L = 0. 6252 

S 2 1 W 21 S 21 u 21 

and increases m 0 to slightly more than unity. Substituting the altered values for weight 
^ 2 

factors into equation (88) yields the new values of as 1310. To this figure must 

be added the increase in variance due to the Russian roulette process as given by equa- 
tion (41): 


A(Ncr^) = p. .w. -S. .r..(l - r. ,)E?. fgl 

v ’ Lj r xj i 3 i] ij' 13 l&j 

M 


- f- i2 w 12^ 1 " r 12^ E 12^1 + ^ 13 W 13^ 1 ' r 13^ E 13^] + ^22 w 22^ “ r 22^ E 22f g ] 

- 38. 9 + 14. 2 + 3610 - 3663 
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II 


Hence, the total variance using almost optimum weight factors is 

Nct t 2 t J = 1310 + 3663 = 4973 

U *l m 2 =l. 004 

which is a figure that is larger than that obtained with straightforward sampling. This 

large value is due to the Russian roulette process practiced on the i=2 j=2 branch which 

2 9 

results in the addition of 3610 to No . The reason the term P22P22 W 22^ " r 22^ z2^ 

2 2 

is much larger than the other terms in A(Ncr ) is that the product p^p^Ef^fg] is much 
larger for the i=2, j=2 branch than for the other two branches in which Russian roulette 
is used. 

2 

In order to eliminate this large value of A(Na ), I/S 22 is taken as equal to 1 giving 
rise to the new values 



s 22 



J- L= _L = 1.875 

S A 

22 u 22 u 22 


and thereby raising the value of 


m 2 to 1. 551. 



m 2 =l. 551 


Now the total variance becomes 


= 1263 


where the contribution to this figure by the Russian roulette process in the splitting step 
is only about 50 in magnitude. 


Case (2) I-I 


Equation (32a) becomes 


Nct h = E ^ij u i u ij T ij - e2 ^ 
i, j 


By table I, for this case, 


1 _ 1 _ _ JiL 

G. E [ T iil’G.. E il 




T ii J 


(91) 


(92) 
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and 


N°n = E 2 [ T ij] “ E 2 [g] = 2194 (93) 

Here the reduction in variance is 

n (4u - = e H] - e2 [v - E ( T ii - E [ T iii)] ‘ 1411 

and the optimum weight factors as given by equations (92) are 


i 

i 

l/“i 

1/G ii 

1 

1 

1. 775 

1. 078 


2 


1. 206 


3 


0. 0682 

2 

1 

0. 668 

0. 640 


2 


1. 744 


3 0 
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E' 


No 


: [ y i <T ii ) ] 


1410 


*1 


(96) 


m., 


m, 


where the correction based on the splitting step of the composite stage is not yet incor- 
porated. 

For = 1, 


/\ 

N * CT uu " 


= E 


m.=l 


2 
T. . 

1 ] 


- E 2 [g] - E 2 [ ri (T..)^ 


= E 


/ \ 2 


/ 

r n \2 


/ \ 2 

(t.. - E.fr. . 1 ) 

V i] lL i] 7 

+ E 

(i<V - E 

[-i (T ii>] ) 

+ E 

[Ejg] - E[g]J 


- 708 + 1304 + 183 = 2195 

The increase in variance due to the splitting step in the 4» stage must be evaluated. 
We get, for m^ = 1, 


i 

1/U, 

l/w t 

l/®i 

1 

0. 1071 

2. 47 

23. 04 

2 

1. 383 

0. 37 

0. 2675 


Again, if Kahn's procedure is followed, then l/s^ should be taken as 23 whereas l/sg 
is left at the value 0. 2675. Applying equation (41), we find 


A(Nct 2 ) = Yj P i w i s i r i ( 1 - r i) E f[g] = P2 W 2 (1 " r 2 )E 2t g ] = 1419 


/\ 

2 


This value, if left unchanged, would more than double the original value Nct^j. Hence, 
we take (I/S 2 ) = 1, which eliminates Russian roulette entirely. Our altered values are 
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i 

j 

1/^ 

l /s i 

1/Wj 


1 

1 

0. 1071 

23 

2. 463 

1. 078 


2 




1. 206 


3 




0. 0682 

2 

1 

1. 383 

1 

1. 383 

0. 640 


2 




1. 744 


3 0 

The value of has been raised to 1.707, and the total variance becomes 


Ncr 


2 


m =1. 707 


E p i w i r i 2<T « ) 


= O^^ 8595 + = 1047 + 98 = 1145 

2.463 1.383 


Case (4) 4>-4> 

The equations are 

2 V - ' 2 2 

Ncr T T = > p.p..w..cr + p.p..w.u..E..[g] 

F i H i] i] g t j 2 -j F i F i] i i] i] L&J 

i, j i,l 

- £ PjW.Effg] + £ P^jEfrgl - E 2 [g] 
i i 

where, by table I, 


1 = E i g 2, 

1 

= E tj [g] , 1 2 hi 

u. E fel ’ 

U- • 

E i^ w., Ei 



1 

1] 

1J 

L g iiJ 



(97) 


( 98 ) 
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and the variance of the sampling process does not depend on the form of the splitting 
factor s^. We get 



E 



m i m 2 


( 99 ) 


As mentioned previously, for optimum results, the splitting step of the first composite 
stage is eliminated so that 


— = 1 , 1 = 1,2 

s i 


(100) 


and m^ is 1. Thus, 4'-4' goes over to I- vp- and equation (99) now becomes 

„2| 




V. 

L g iJ 


m r 


1 133 

m„ 


( 101 ) 


For m„ = 1, 


/\ 

, T | 2 2 

N|CT uu - a I* 



i = E 

r 2 1 

T. . 

m 2 =1 y 


_ 1 ] 


- E 2 [g] - E 2 


(7 


g ij. 


= E 


/ 

r 1 

\ 2 1 


/ \ 

■ E 

r 

crq 

M* 

«— 1. 

1 

) 

+ E 

- E fg]) 


= 1673 + 799 = 2472 


and 
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i 

j 

1 

1 

1 



u.u. . 
1 lj 

w. . 
i] 

"ij 

1 

1 

0 + 

2. 97 

OO 


2 

4. 32 

1. 485 

0. 3435 


3 

0. 2593 

0. 0594 

0. 229 

2 

1 

0. 432 

0. 594 

1. 372 


2 

2. 593 

0. 297 

0. 1145 


3 

0 

0 

0 


Again, by Kahn's procedure, we take l/sgj equal to unity. In addition, to avoid 
an unduly large variance due to the splitting method, I/S 22 is also se t equal to unity. 
Hence, l/w 2 j = 0- 432 and = 2. 593 and these values correspond to a value for 

m 2 of 1. 781. 

The formula for total variance now becomes 



m 2 =l. 781 


Z 1 &..W ..( 7 ^ +\ j©..w..s..r..(l - r..)E?. fgl 

r iJ 13 g ij 13 13 13 13 13 13 L J 

i,j 


= 1063 + 38 = 1101 


SUMMARY AND CONCLUDING REMARKS 

General formulas have been developed for the variance reduction obtained under 
different sampling schemes for the situation where the two techniques of importance 
sampling and splitting and Russian roulette are employed. Optimum biasing proce- 
dures have been determined for any number of independent random variables and the 
results are shown in tables I and II. 

It is found that, when the random variable under consideration is non-negative, the 
minimum variance with the least number of members in the sample is obtained if the last 
random variable is sampled in a combination importance-splitting sampling manner and 
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I 


all preceding variables are sampled employing adjoint biasing. A short example has 
illustrated the possible increase in variance that may occur if indiscriminate Russian 
roulette practices are used. 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, September 22, 1970, 

129-01. 
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APPENDIX A 


A SPLITTING TECHNIQUE 


We are to consider the splitting method in which, for each of the original 
X +1 ' 


n. 


members, < j 1 * > members are generated, the last member having a weight equal to 

lij J 


X 


1 L that of the first i j 1 > members. Thus, an effective splitting factor equal to 


U 


n _i +r 

7 l ' 1 


l] 


-L = I..+r.. 
e U U 

S ij 


>is obtained. 


The random variable to be considered is that given in equations (44) and (45) of the 
text and repeated here: 


£ 

ij 


Y => Y ii 


(44) 


where 


+ r. . 
i] 


[\h 

v. • 

x 2h 

Z ( g «Ki +r ‘ 

Z 

["lr 1 

a 2r 1 

\h 

V . . 

X 2 J 2 

Z ( g «)“u +r ‘ r 'J 

z 

01 12 _1 

a 22~ 1 


la 


21 


l a 


22 


(45) 


The mechanics of carrying out the sampling procedure are as follows. A sample of 
size N is first picked from the x. -population in accordance with the probability distri- 

1 4-U 


,th 


bution function p^/u^ and n. designates the number of members possessing the i 
value set of x.. These n. members generate u, = I.n. members, each having an 

X 1 1 ^ 11 

associated weight of w^, and u i = n^ members, each having a weight of w-r^. The 
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sample members are further subdivided in accordance with their x 0 -values where the 

probability distribution function is p../u.,. Thus, of the members possessing the i Ln 

th *-3 

value set of x 1 and the j value set of x 0 , there are n. • members with a weight of 

1 ^ 1 jJ 

w.u. where \ n. - = v. and n. , members with a weight of w-r.u. where 
1 1 ijJ ^3 111 

\ n. • = v. . Finally, each of these members in turn is split into (I. +1) members 
L-i x 2 3 l 2 
j 

where the last member has a weight r^ that of the others. Hence, for each i and j, 

there are four groups: (1) v. ■ = I..n. , members with weight w.. , (2) v. ■ = I n 

1 1 3 1 13 x l 3 13 1 2 ] 1 iJ l 2 3 

members with weight w.,r., (3) v. ■ _ n. • members with weight w..r.., and (4) v. ■ 

i] 1 1^2 “ ^3 i] i] J 2 ] 2 

= n. . members with weight w, r,r... 

We wish to determine the variance of Y. Proceeding as in equation (28), we have 


2 

ct Y = 


Z e M + Z 


E Y..Y.., 
1 i] 13’ 


l » 3 


b3,3' 

(j/j') 


z 

i,j,i',3' 

dA') 


Y..Y.,., 

13 13 


E 2 [g] 


(Al) 


The terms of the first summation on the right side of this summation may be written as 


Y 2 

ij 


w 2 

Je 

N 2 


li ! ( Si j) 1 1 +r i y i 2 j 1 ( g ij)21 +r ij* , i 1 j 2 ( g ii)l2 + r i r ij ^i 2 3 2 ( g ij) 22 


(A2) 


There are ten terms to be evaluated, four squared terms and six cross products, and 
each can be handled by the methods employed previously. For example, 


E 


ViN 


2 

11 


= E 


Si 


EfJg] 


\h\ -j j 131 


E 

g ij 


y ii3 


1 J 1 


+ E^[g]E 


V * 


l l 3 l 


9 Pi Pii 9 9 

= 0 N — I. I.. + E .[gll . E 
g.. 1 11 n | ] n 

13 u i u ij 


' 2 
V 


57 


L 



and 


E 2 Vi 1 j 1 fo)ll , 'l 1 i 2 ( g lj)l2 -»>-l|*SW“[-l 1 j 1 *’l A 


= 2r.,E 2 fg]I..E n? , 
i] i] L6J i] 


where 



E nr • = N‘ 


2 Pi Pj / Pi \ 2 p ii p i p ii 1 p ii 
1 — + N — 1 - -i ) I -ii + N-i I. _U 1 - _il 


2 u- \ u. / 1 2 u. 1 u. . \ u. • 

i \ 1/ u A - 1 1] \ 1] 


Evaluating each of the terms entering into equation (A2), we find that 


E k] = E f<?ij - i {> - 1 1 - s i r i (i - r i>] t 1 - s ii r ii (i - r ij>]} 


-^Wi s i u ii( 1 -^ r i (1 ' r i )E ii [E) 


where the expression for E r.j is given in equation (29). Continuing with the computa- 
tions, we find that, for j/j', 


Vij . W 'V :; E ( 1 'i 1 i 1 ( E ii)n + r i‘'i 2 j 1 («i]>21 + Vi li 2 (Sii)l2 + r i r ij ^2^)22 


'ijji( g ii') ii + r i l 'i 2 ii( g ii l ) 2i + r tr V 2 ( g ir ) 12 + r i r iJ^ 2 J2( ei r ) 22 


• E i’ij 3iy + 5 p i p ij p ir w i s i r i (1 - r i) E ij[sJ E ij'[s] 


and, for i/i', 
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E [ Y ii Y i'i'] = E 

where the expressions for E and E are S* ven i n equations (30) 

and (31), respectively. 

Substituting equations (A3), (A4), and (A5) into equation (Al) results in 


Nct 


2 

Y 



Z p.p..w..a 2 f 1 - \l - s.r.(l - r.)l |~1 - s..r..(l 
i] i] gy l L i i v rj L i] i] 

M 



I 


bl 


p ii 

p i w i s i r i ( 1 - *.) — 

U ij 



(47) 


Similarly, it can be shown that, for the splitting method considered in this appendix, 
if there were only one sampling stage, then the change in variance would be given by 

A(Nct 2 ) = - p i w i CT g. s i r i( 1 ' r i) ( A6 ) 

i 

and, if there were three stages of sampling, then 


A(Nct ) = 


I 

i,j,k 


p i p ij p ijk w ijk CT gi . k l 


{l- [ 1 - s i r i< 1 ‘ r i)][ 1 


s.r..(l 
i] i] v 


V. 




Z p.p..w..s..r..(l - r..) (u. E... [g] - E..[gl) 

H i F i] i] i] i] v ij 7 l i]k i]k L6J i] L&J / 

ijk x 

i,l,k 

- ^ PiWiVid - r.) (u^E-.fg] - E.[g]^ 

i,i 


(A7) 


The formula for any arbitrary number of stages should now be readily apparent. 
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APPENDIX B 


PROOF OF INEQUALITY 

We wish to show that 

E a > E a 
g ij [ g ijk 

We can write 

(J g i _ ( g ijkl " E i^) P ij P ijk. . . 

3> k, 1, . . . 

'S - ' l 2 

= Z, [fai*. . . - V gl ) + (v g J- E i [g] ) p ij p ijk. . . 

j»k,l, . . • 

X' : 'l X ( eiik1 ' • ' ' E » M ) p lik p ukl. . . 

i k, 1, . . . 

+ X P «( Ei i [g| ' E ‘ [g| ) 2 

j 

= X Pii f + *' Eii[gl ’ ElW )j 

= E i % + ( E ij [gl - E il g l) 2 

and similarly 

ctJ =E.. cr 2 + (E... fgl - E..[gl \ 2 

g ij 1] L g ijk \ 1]k 13 / 

Hence, 


(61) 


(Bla) 


(Bib) 


60 




(B2) 


which can alternatively be shown by proceeding as follows: 


E. 


g ijk 




W* g ijk 



Pi i P ijk a r?ijk 


2 


^^ p ij p ijk cr g i j k p ij tP ij'k' 
j,k j’,k’ 


y', p ij p ij' p iik p ij'k’ CT g ijk CT g irk , 

j,k, j’,k' 


z 

j> k, j’,k' 


p ij p ij ,p ijk p ij'k' 


,cr -a o 

L g ijk g ijk g ij'k'J 


■i z 


P ij P ij ,P ijk p i3k’ 


j,k, j\k' 


a - a 
|_ ^ijk Sij' k > 
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Thus, 




1/2 


> E. 


cr 


T g ijkJ 


(B3) 


Substituting inequality (B3) into inequality (B2) yields 


cr 

> E 

E. 

a 


= E 

a 

L g d 


l 

_ g ijk_ 



. g ijk 


(B4) 


and inequality (61) has been proved. It may be remarked that inequality (61) merely 
demonstrates the obvious fact that measurement of an additional two correlated variable 
sets in general reduces the expected standard deviation. 
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APPENDIX C 


VARIANCE RELATION FOR NONOPTIMUM SAMPLING 

In this appendix, we will determine, for a 'k-I-S-I sampling, the optimum values 
of N, m^N, and m^m^N that minimize the variance subject to the constraint given by 
equation (74). It will be assumed that the weight factors m, w^, u^, w^.^, and u.^ 
are not at the optimum values as given by equations (72). 

Using equation (49b) where k = 4 and the relations 


w. . = u. .w. 
U i] i 


u. =1 
Uk 




w ijkl = u ijkl w ijk J 


(Cl) 


that hold for a 'T> - 1 - S - 1 sampling, we obtain 


tvt 2 

Ncr = 


E 

i> i) K i 


/S ijkl w ijk u ijkl T ijkl 


- E / 2 ijk' , iik E , 2 jk[ g l + E ^iik w i u ij E i 2 jkfe] 

i,j,k i, j, k 


yy /»iWi E i 2 f g ] + ^^E^g] - E 2 [ g ] 


(C2) 


We can write each of the weight factors as a sum of its optimum value and its deviation 
from optimum 


u. = u. + Au. 

li l 


w^ = w^ + Aw. 


w iik = "ijk + A »ijk 
•ijki ,: iiki • a 


(C 3) 


63 


where the optimum values are given in equations (72). Substituting equations (C3) into 
equation (C2) gives 


I 


Na = > /* ijkl 


w... + Aw.., Au..,, + u..., Aw.., 
i]k i]k I ljkl i]kl 13k 


T. 


ijkl 


ij, k i, j, k 

- ^ ^ i (Aw i )E i 2 [g] + ^ i (Au i )E i 2 [g] 


w i + Aw iy Au ij + u ij Aw i 


E i 2 jk f8j 


E 2 

r,(e.. 


1 ij 


E 


m. 


r ijk (T ijkl ) 


m l m 3 


(C4) 


2 

The sum of the last two terms on the right side of equation (C4) is Nct (see eq. (73)) 

2 

and, hence, in the event that all deviations from optimum go to zero, Nct reduces to 
2 

Nct as it should. /s 

2 2 
We would like to show that Nct is the smallest value that can be obtained for Ncr 

by establishing that the sum of the terms on the right side other than the last two terms 

2 2 

is larger than or equal to zero. In order to accomplish this, let us define A , B , 
and C 2 by 


A' 


-z 

U,k, 1 


^ijklKjk Au ijkl +U ijkl Aw ijk T ijkl 


-Z'« 


i3k( A ^k> E i3kte] (C5) 


ij, k 


B 


= Z au ii + “ij aw i) E ijk[s] - Z ^i (aw i )E i 


[g] 


i, j,k 


(C6) 


C 2 = 


Z^i <aU i )E i 


[g] 


(C7) 


so that 


XT 2 .2 „2 „2 

Ncr -A +B + C + 


y i^ e ij) E [ y ijk( T ijkp] 


m. 


m l m 3 


(C8) 
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1 


With the use of the constraining conditions on the weight factors as given by equa- 

2 2 2 

tions (71), it can be shown that each of the quantities A , B , and C is larger than or 
equal to zero. 

o 

Let us treat the quantity C in detail to illustrate the procedure. By equations (71), 
satisfies the equation 



i 


and, consequently, 



i 


(C9) 


(CIO) 


Equation (C9) can be written 


where 




Introducing equations (CIO) and (C12) into equation (Cll) gives 



i 


(Cll) 


(C12) 


(C13) 


which is the constraining condition that (Au^ must satisfy. Now equation (C7) can be 
expressed as 


C 2 = E 2 [g] 



Pi A Ui 


~2 

u i 


(C 14) 
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II 


where equation (72a) has been employed. We are now in the position where we can mod- 
ify equation (C14) by using equation (C13) to yield 


C 


2 



= E 2 [g]J~^-hj (iu.) 2 = E 2 [g] (Au.) : 


(C15) 


Equation (C15) is the desired relation that shows C > 0. 

2 2 

Treating the quantities A and B in like manner we have 


a2 = Y ^iikl w ijk< au i)kl> T «M + Y 




i,j,k 


p ijkl u ijkl r ijkl ” E ijk^ 


\ P ijkl w ijk^ Au ijkP ~ 2 + /Z ^ijk^ Aw ijk^ y ijk^ T ijkl^ 

/ U--i,i j ; l. 



ijkl i, j,k 


i,j,k, 1 


z 

i> j> k 


t 1 ijk w ijk E ijZ T ijkl^ 



' p iikl p 


U ijkl 



ijkl 


(Au ijkl> 


u ijkl u ijkl 


1 



^ijk 


- <Aw ijk> 


w ijk w ijk 


I 

i, 1, k, 1 


^ Mp »ijk T f i ki< A “i i id) 2 + y ^-^fjk^ijkix^ijk' 2 

U ijkl L^i w ijk 

i, j,k 


(C16) 
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and 



We see that minimization of the variance as given by equation (C8) with respect to 
the variables N, m^N, and m^m^N subject to the constraint expressed by equation (74) 
no longer yields zero as the optimum value of N but instead gives 
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N 


T( A 2 + B 2 + C ^ /^i 


\ 1/2 

a^A 2 + B 2 + C 2 j + a 2 E|r i (e ii )| + agE| 


: [ y i (6 U>] +a 3 E [njk^iikl)] 


(C18a) 


The optimum values of m-jN and m^m^N become 

TE 

1/2 


mjN ] = — 


a, IA 2 + B 2 + C 2 ) + a 2 E[ r .( £ij )] +a 3 E[y ijk (T ljkl )_ 


(C18b) 


m^mgN) = 


,2 . „2 . ^ 


TE 


1/2 


[ y ijk^ T ijkl^_ 


a,lA- + B- + C-i + a 2 E ^(e..)] + agE [y^r^) 


(C18c) 


and the minimum variance is 


2 \ 1 

a I . = — 

/min T 


a. 1 [A 2 + B 2 + C 2 ) +a 2 E 


2 

!/2 .. 

) +a 2 E [ y i (e ij)J +a 3 E [ y ijk (T ijkl ) J 


(C 19) 


o \ / / 2 \ 

It is to be noted that ^cr } - n /because yr y m i n * s n °t minimized with respect to 

the weight factors. To be redundant, if the weight factors take on their optimum values, 

then A 2 + B 2 + C 2 = 0 and (a 2 ) min = (a^j. 
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APPENDIX D 


A 2 

2 2 2 
a l ,a 2’ a 3 

B 2 

,2 k 2 .2 

b 1 ,b 2’ b 3 

C 2 

c ^l’ c I2’ c S3’ C I4’ C g 

E t □] 
e[D/0] 

EiO 

VO 1 

f(EE) 

f(D/0) 

g 


L 

m i 

N 


SYMBOLS 

quantity defined by eq. (C5) 
coefficients appearing in eq. (74) 

quantity defined by eq. (C6) 
coefficients appearing in eq. (77) 

quantity defined by eq. (Cl) 

coefficients defined immediately before eq. (74) 
expectation value of EE 

expectation value of EH given that O is fixed at a 
certain value 

j_j_ 

expectation value of EE] given that Xj takes on its i n 
value set 

expectation value of EE given that x. and x 0 take on 

J.L J.|_ X d 

their i Ln and j value set, respectively 
probability density function of I I 

probability density function of EE given that O is fixed 
at a certain value 

function dependent on a number of sets of random variables 

integral portion of the reciprocal of s^s^,. . re- 
spectively 

pure importance sampling stage 

importance sampling step in a composite sampling stage 

quantity employed when using the Langrangian- multiplier 
technique of minimization subject to constraints 

th 

magnification factor for i stage 

original number of members in sample 

number of sets of variables upon which g is dependent 
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n. 

1 


n. . 
i] 


P ii 




m 


Pi 


13 


P ijk 


r., r. 

1’ i] 


S 

S T 

s(x 1 ,x 2 , 


s. 

1 


S ij 

T 

F 

U 

u(x p x 2 , 


u. 

1 


u. . 
1] 


w(XpX 2 , . . . ) 


number of members in sample possessing i^ value set of 
Xj immediately after importance sampling on x^-variable 

number of members in sample possessing i^ 1 value set of 
til 

Xj and j value set of x 2 , immediately after impor- 
tance sampling on Xg- variable 

probability that takes on its i^ 1 value set 

conditional probability that x 9 takes on its j value set 

u j 1 

given that x^ has taken on its 1 value set 

conditional probability that Xg takes on its k^ value set 
given that x^ and x 2 have taken on their i^ 1 and 
value sets, respectively 

J.L rr 

probability that x^ and x 2 take on their 1 and j 
value sets, respectively 

probability that x^, x 2 , and Xg take on their i^ 1 , j^, 
and k^ 1 value sets, respectively 


fractional portion of reciprocal of 
spectively 


3 i’ S ij ’ 


re- 


pure splitting sampling stage 


splitting step of a composite sampling stage 


splitting factor dependent on variable sets Xj,x 2 , . . . 

splitting factor for x^- stage of sampling 

splitting factor for Xg- stage of sampling 

constraining quantity defined by eq. (74) 

binomial variable occurring in analysis of Russian roulette 


unit sampling stage 


importance sampling weight factor dependent on variable 
sets Xj,x 2 , . . . 

importance sampling weight factor for x^- stage of 
sampling 

importance sampling weight factor for Xg-stage of 
sampling 

overall weight factor dependent on variable sets 
x 1 ; x 2 ,. . . 


70 



w. 

1 

overall weight factor for x^-stage of sampling 

y 

overall weight factor for x^- and Xg-stages of sampling 

X 

collection of sets of random variables upon which g is 
dependent 

x i 

i^ set of random variables 

Y 

random variable defined by eq. (44) 

Y U 

random variable defined by eq. (45) 

Z 

random variable defined by eq. (16) 

Z ii 

random variable defined by eq. (17) 

J 

random variable defined by eq. (24) 

3 ij 

random variable defined by eq. (25) 

a . 

l 

absolute value of E.[g] 

a .. 
i] 

absolute value of E.j[g] 

y 2 0),y i 2 (n),y 1 2 j (0.- ■ 

quantities defined by eqs. (67) 

A(Q) 

change in | 1 

6 2 , 6 2 ,5 2 .,. . . 
’ i’ i]’ 

quantities defined by eqs. (66) 

2 2 2 
£ ’W ‘ ‘ 

quantities defined by eqs. (54) and (65) 

K 

number of stages 

A 2 

Lagrangian multiplier 

V 

number of members in sample immediately after splitting 

2 

a 

variance 

2 

°n 

variance of Q 

°g/ x l’ x 2 

variance of g when x^ and Xg are fixed at certain values 

(A ■ 

\ / min 

minimum variance appearing in appendix C 

2 2 2 

quantities defined by eqs. (64) 

* 

composite sampling stage 

*23 

composite sampling stage for variable sets x^ and X 2 
simultaneously 
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□,o 

Subscripts : 
i 

i 

k 

a,p 

Superscripts: 
' (prime) 


any quantity 

Xj fixed at its i^ 1 value set 

X. L 

Xg fixed at its j value set 

j-i_ 

Xg fixed at its k value set 
dummy indices 

differentiates from unprimed quantity 
optimum value 
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TABLE I. - OPTIMUM RECIPROCAL WEIGHT FACTORS (g NON-NEGATIVE) 



E [g] E.[g] 


TABLE II. - OPTIMUM RE CIPRICAL WEIGHT FACTORS (NO' RESTRICTIONS ON g) 



NASA -Langley, 1971 
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